This effort is also part of the VAIA FRONT project - FRom lessong learned to future Options .

Summary

We have to assign land-cover attributes to nodes over a large area for simulating the impact of strong winds over forest areas.

Several layers are processed to aggregate information over nodes of a regular grid.

Data were processed with Google Earth Engine (GEE) and R.

Code and outputs

Code is reported below (snippets) and available in Github .

Outputs are available at Github page in the folder “out”.

Forest DAMAGES layer

Training and validation was done over areas that were defined over areas that were damaged from the VAIA storm event.

The data is extracted from the windthrow database that was created for the paper by (Forzieri et al. 2020) linking to an open database containing polygons of damaged areas. It must be noted that some damaged areas were left out because only slightly damaged, or with small areas or outright missed during the visual interpretation process. Therefore we should be careful to conclude that all areas not covered by the polygons are un-damaged. Nevertheless these data provide us with samples of damaged areas.

There are 7416 polygons in the dataset: 5235 have information on the degree of damage, expressed between 0 and 1, 2181 have NULL value (expressed as -9999).

Aggregating data at nodes

A regular grid with 1 km spacing was used over the study area, for a total of 97524 nodes.

N.B. A nice alternative could be to have an hexagonal grid, which has several geometric advantages.
Cannot convert our regular grid to hexagonal grid because hexagonal requires offset at alternate rows.

Grid is obviously not aligned to a regular latitude longitude grid as distance between nodes would decrease moving away from the equator.

A custom projection Lambert Conformal Conical projection was designed (Marika Koukoula) with the following PROJ

myproj<-"+proj=lcc +lat_1=45.827  +lat_2=45.827  +lat_0=45.827 +lon_0=11.625 +x_0=4000000 +y_0=2800000 +ellps=GRS80"

Each node was expanded to a square

  st_bbox_by_feature = function(x) {
    x = st_geometry(x)
    f <- function(y) st_as_sfc(st_bbox(y,), crs=myproj)
    do.call("c", lapply(x, f))
  }
  
    
  nodes.myproj.buffered<-st_buffer(nodes.myproj, 500, nQuadSegs = 1)
  squares<-st_bbox_by_feature(nodes.myproj2)
  squares <- squares  %>% st_set_crs( myproj)

The square lattice was used for further processing with Google Earth Engine and with R.

Google Earth Engine processing

The 97524 square polygons were uploaded to Google Earth Engine for processing raster products derived from satellite data.

The grid was used to extract of the following two files:

The following code was used in GEE.

var getCentroid = function(feature){
  return feature.centroid();
};  
 
var slp = ee.Terrain.slope(srtm).clip(domain_buffered.geometry().bounds()) ;
var asp = ee.Terrain.aspect(srtm).clip(domain_buffered.geometry().bounds()) ;
var dem =  srtm.clip(domain_buffered.geometry().bounds()) ;
 
var image = hansen.select(['treecover2000']).rename(['cCov'])
            .addBands(canopyHeight.rename(['cHgt']))
            .addBands(slp.rename(['slp']))
            .addBands(asp.rename(['asp']))
            .addBands(dem.rename(['dem'])) ;


var  treeCover =  image.reduceRegions({
  reducer: ee.Reducer.percentile([10,25,50,75,90]),  
  collection:domain_buffered 
}); 
 

var  unBuffer = ee.FeatureCollection(treeCover.map(getCentroid));
 
 Export.table.toDrive({
    collection: unBuffer,
    description:'mapRedVars',
    fileFormat: 'SHP'
  });

R-CRAN processing

Processing 97524 squares intersecting them with 121924 polygons of forest category and 7416 polygons of damage areas too 120 seconds with 14 CPU threads parallel tasking (not sure if this is the correct terminology - not a parallel guru).

Variables

Table with variables: asp=aspect, slp=slope, dem=DEM, cHgt=canopy height, cCv=canopy cover, mainCat=main category, mnCt_Cv=main category cover, dmg_Cov=damage cover, dmg_wCov=damage weighted cover

All cover values are in m^2 and can be converted to relative cover dividing by the square area (1 km^2).

Column Name Type Min Mean Max
asp_p10 numeric 0 71.05 336.79
asp_p25 numeric 0 116.13 341.25
asp_p50 numeric 0 174.04 351.06
asp_p75 numeric 0 232.57 355.19
asp_p90 numeric 0 273.54 358.60
cCv_p10 numeric 0 10.84 95.00
cCv_p25 numeric 0 22.00 100.00
cCv_p50 numeric 0 36.69 100.00
cCv_p75 numeric 0 49.42 100.00
cCv_p90 numeric 0 58.46 100.00
cHgt_10 numeric 0 9.95 36.00
cHgt_25 numeric 0 12.13 36.00
cHgt_50 numeric 0 15.35 40.00
cHgt_75 numeric 0 18.31 40.00
cHgt_90 numeric 0 20.04 40.00
cat_Cov integer 0 597 686.44 2 999 999.00
dem_p10 numeric -11 1 071.91 3 566.37
dem_p25 numeric -7 1 118.29 3 611.00
dem_p50 numeric -6 1 183.86 3 701.17
dem_p75 numeric -6 1 253.49 3 794.00
dem_p90 numeric -6 1 308.16 3 919.00
dmg_Cov integer 0 90 426.05 1 074 988.00
dmg_wCv integer 0 66 670.52 963 963.00
mnCt_Cv integer 0 427 530.22 2 820 898.00
slp_p10 numeric 0 10.54 41.61
slp_p25 numeric 0 14.41 48.66
slp_p50 numeric 0 19.01 65.91
slp_p75 numeric 0 23.60 73.24
slp_p90 numeric 0 27.48 82.84

Map of subset of study area

Forest area and forest categories

Italy is divided into regions and each region into provinces. We have a full cover of forest categories of two regions, Veneto and Friuli Venezia Giulia, and the province of Trento. Categories slightly differ in classification, but discrepancies were merged into a single catogory.
Classification can provide us with a binary result: positives/negatives = damaged/undamaged forest areas - with the following alternative outputs per each unit:

  • true positives = classified as damaged correctly
  • false positives = classified as damaged but not true
  • true negatives = classified as undamaged correctly
  • false negatives = classified as undamaged but are actually damaged

Forzieri, Giovanni, Matteo Pecchi, Marco Girardello, Achille Mauri, Marcus Klaus, Christo Nikolov, Marius Rüetschi, et al. 2020. “A spatially explicit database of wind disturbances in European forests over the period 2000 - 2018.” Earth System Science Data 12 (1): 257–76. https://doi.org/10.5194/essd-12-257-2020.